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ABSTRACT 


As  linear  programming  techniques  find  applications  in  more 
diverse  fields,  the  problem  of  solution  time  becomes  increasingly 
important.   A  variation  of  the  revised  simplex  algorithm,  in  which 
the  constraints  are   added  in  a  step-by-step  fashion,  is  investigated 
as  a  potentially  faster  solution  technique.   A  computational  pro- 
cedure, coded  for  the  IBM  3&0  computer,  is  developed  to  compare  this 
algorithm  with  the  standard  two-phase  revised  simplex  algorithm. 
A  limited  number  of  problems,  including  several  randomly  generated 
problems,  is  solved  by  each  of  the  two  methods.  The  resulting  com- 
parison of  solution  times  indicates  that  a  significant  improvement 
is  obtained  by  the  use  of  the  procedure  of  step-by-step  addition  of 
constraints. 
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1 .   Introduction 

Increased  utilization  of  linear  programming  techniques  has  led  to 
the  formulation  of  problems  of  sufficient  size  to  tax  the  storage 
capabilities  of  many  computers  now  in  use.  The  problem  of  storage 
capacity  may  be  alleviated  by  auxiliary  storage,  with  the  accompanied 
requirement  of  access  delays  which  increase  solution  time.   Even  with 
sufficient  storage  and  a  coded  procedure  large  enough  to  accommodate 
a  problem,  the  computation  time  may  be  such  that  solution  costs  ap- 
proach a  budget  limitation.   For  example,  the  CEIR  LP/90  program  for 
the  IBM  7090  computer,  which  can  accommodate  512  constraints  and  an 
unlimited  number  of  variables,  or  the  Philco  Corporation  LP  2000 
System  for  the  S-2000  computer,  which  can  accommodate  2500  constraints 
and  an  unlimited  number  of  variables,  might  take  several  hours  to 
solve  a  large  problem  at  a  cost  of  several  thousand  dollars. 

Improvement  in  the  solution  time  might  be  accomplished  by  im- 
provement of  the  hardware  or  by  an  improvement  in  the  mathematical 
procedure.  The  hardware  improvements,  which  are  becoming  available 
in  the  newer  computers  such  as  the  IBM  3&0,  are   increased  computa- 
tion speeds  and  faster  access  to  external  storage. 

Solution  procedure  improvement  has  been  attempted  by  several 
methods.   Some  attempts  are  being  made  to  take  advantage  of  special 
structuring  of  the  problem.  These  techniques  are  directed  toward 
partitioning  or  decomposing  the  problems  into  manageable  sub-problems. 
Another  approach  to  improved  solution  time  is  the  Primal-Dual 
algorithm  which  was  developed  by  Dantzig.  [3] 

This  project  has  been  directed  towards  another  possible  method 
for  improving  solution  procedures.  An  appraisal  is  made  of  the  step- 


by-step  addition  of  constraints  to  the  two-phase  revised  simplex  pro- 
cedure as  a  possible  solution  technique  which  would  serve  to  improve 
solution  times.   This  method,  called  here  the  "step-by-step  addition 
of  constraints  (SSAC),"  exploits  the  advantage  of  obtaining  a  rapid 
solution  to  a  small  sub-problem  and  then  moving  from  one  optimal 
solution  to  an  optimal  solution  of  a  slightly  larger  sub-problem  as 
the  remaining  constraint  equations  are  added  one  at  a  time.   In  the 
earlier  stages  of  the  solution  procedure,  the  size  of  the  matrices 
used  in  the  iteration  procedure  would  be  relatively  small,  compared 
to  the  complete  problem.   By  exploiting  the  advantage  of  multiplying 
smaller-dimensioned  matrices,  with  the  attendant  shorter  computa- 
tional times,  a  shorter  overall  solution  time  might  be  achieved  in 
spite  of  the  fact  that  a  greater  number  of  iterations  would  be 
required. 


2.   Notation 

The  notation  used  throughout  this  paper  was  chosen  to  correspond 
to  the  notation  most  frequently  used  in  linear  programming  texts. 

(1)  Upper  case  letters  represent  matrices. 

(2)  Lower  case  letters  represent  column  or  row  vectors. 

(3)  Subscripted  lower  case  letters  represent  elements  of  row  or 
column  vectors. 

{k)      Tableau  notation,  as  illustrated  below,  is  similar  to  that 
used  by  Dantzig.  h] 

Z    W    Xr   X^   X2   Xo   x^    b 


*-«- 


-2    0  -k        0 


.3  -k     „6  -1  -2k 


0 


3  k     ©  1  I  2k 


*  -  2k/e  =  k 


'Rl 


The  columns  corresponding  to  basis  vectors  are  indicated  by  a  dot. 
The  pivot  column  (or  row  in  the  case  of  the  dual  simplex  algorithm) 
is  indicated  by  an  arrow.   The  pivot  element,  as  determined  by  the 
appropriate  minimum -0- cr  iter  ion,  is  indicated  by  circling  the  element, 
and  the  basis  variable  which  is  then  to  be  driven  out  is  indicated  by 
circling  the  associated  dot. 

The  tableau  rows  are  referred  to  as  the  "z",  "w"  or  Ri"  row  where 
i  is  an  integer  corresponding  to  the  sequence  in  which  the  constraints 
appear.   That  is,  MR1"  refers  to  the  first  constraint,  "R2"  refers  to 
the  second  constraint  and  so  on. 

Specific  notation  which  will  receive  repeated  use  in  this  paper 
inc ludes: 

m   number  of  constraint  equations 


n   number  of  variables 

A   m  x  n  matrix  of  coefficients  of  the  constraint  equations, 

having  elements  a-° 
1  J 

B  m  x  m  matrix  of  basis  vectors 

B  inverse  of  the  basis 

P-  the  m-dimensiona 1  vectors  which  make  up  the  B  matrix 

x  an  n-dimensional  column  vector  having  elements  x^ 

c  the  n-element  row  vector  of  cost  coefficients  having  elements  c,- 

b  the  m-element  column  vector  of  the  right-hand  side  of  the  con- 
straint equations  (requirements  vector),  having  elements  b- 


3.   Formulation  of  the  Problem 

The  general  linear  programming  problem  is  stated  as: 
Maximize 

ex 
subject  to: 

Ax  =  b, 
and, 

x  >  0. 
For  the  algorithm  to  be  investigated,  we  define 

z  =  ex. 
We  can  then  rewrite  the  problem  as: 
Maximize 

z, 

subject  to: 

Ax  =  b, 
z  -  ex  =  0, 
and, 

x  >  0. 
To  obtain  an  identity  matrix  to  begin  the  two-phase  revised 
simplex  procedure  we  add  artificial  variables  xn+],  xn+2>  ....,  xn+m 
to  the  constraint  equation  with  xn+,-  ^  0  for  j  =  1,...,m. 
I f  we  def  i  ne 

w  -  -xn+]  -x-n+2  "•••~xn+m 
and  then  maximize  w,  we  will  drive  the  artificial  variables  out  of 
the  basis  and  either  obtain  an  initial  basic  feasible  solution  or  an 
indication  that  the  problem  is  infeasible.   The  process  associated  with 
maximizing  w  is  usually  referred  to  as  Phase  I. 
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During  the  second  phase  (Phase  II),  the  problem  is  stated  in  the 
f ol lowi  ng  form: 
Maximize 

z 
subject  to: 

a11x1  +  a12x2  +  +  a1nxn  +  xn+1  =  b1 

a21x1  +  a22x2  +  +  a2nxn       +  xn+2  =  b2 

am1x1  +  am2x2  +  +  amnxn  +xn+m  =  bm 

"c1x1  "  C2X2   cnxn  +  z  =  ° 

xn+1+xn+2+--+xn+m+  w  =  ° 

and, 

Xj  >  0   for  j  =  1,  2,  ...,  m. 
Because  of  the  w  equation  we  do  not  have  an  identity  matrix  to 
use  as  a  starting  basis.  Therefore,  we  subtract  from  that  equation 
each  of  the  other  equations,  as  appropriate,  to  remove  the  xn+j,  ..., 
xn+m  vafiables.  We  obtain  the  following  form  for  the  w  equation: 

w  +  am+2,1x1  +  am+2,2x2  +  •••  +  am+2,nxn  *  bm+2 
where, 

m 


am+ 

i  =  1 


and, 


m 

b*2      =  -  E  bi 

l  =  1 
for   those   artificial    variables,    i,   having   non-zero  prices.      Our   problem 
is    now   in  a   form  such   that   we   have  an   identity   matrix   to   use  as   our 
starting   basis. 
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Maximize 


subject    to: 


3  1  iX  i     "T*     8  i  oXa     *•" 


3o  i  X  i     ">"     3^aX<j     "t* 


+   alrx     +   x    .  , 
I  n   n  n+ I 


+  a~x_  +   x 

zn   n 


n+2 


=  b 


=  b. 


am1x1  +  am2x2  + 
-c]x1  -  c2  x2  - 


+  a  x 
mn  n 


"  cn  xn 


+  x 


n+m 


+  z 


am+2,1xl+am+2,2x2  +  '"  +  am+2,nxn 


and, 


+  w 


=  b 

i=  0 
=  b 


m 


m+2, 


x  >  0  for  j  =  1,2,  ...,  n+m. 
If  we  let  a^.i  •  ==  -c  •  ,  the  problem  resolves  into  the  three 
matrices  used  for  computations  in  the  revised  simplex  procedure: 


A  = 


a11  a12 


a21  a22 


aml  am2 


1n 


'2n 


mn 


m+ 1 , 1 m+ 1 ,  n 

am+2,1 am+2,n 


1  0 
0  1 

• 

0 


0 


,  an  identity  matrix,  and 


b  = 


m 


0 

bm+2 
For  simplification  in  handling  the  matrices  when  coded  in  FORTRAN, 

they  will  be  re-arranged  so  that  the  objective  function  and  the 

modified  w  equation  appear  in  the  first  two  rows  rather  than  the  last 

two: 


m+1,1 


W,1 


11 


'21 


m+1  ,  n 


m+2,n 


In 


'2n 


'ml 


mn 


3m+2 
B"   =  I,  and   b  = 

L 
m 

In  the  revised  simplex  procedure,  the  original  A  matrix  and  b  vector 
of  the  full  starting  tableau  ("original"  tableau)  are   used  at  each 
iteration  with  the  inverse  of  the  current  basis  matrix  to  determine 
certain  unknown  elements  of  the  current  tableau. 


)k 


When  a  constraint  is  added,  the  original  A  matrix  and  b  vector  are 
changed.   The  size  of  each  is  increased  by  one  row,  and  the  A  matrix  is 
increased  by  one  column.  The  basis  matrix  is  therefore  increased  in 
size  by  the  addition  of  a  row  vector,  Q    ,  corresponding  to  the  elements 
of  the  new  constraint  which  are  in  the  current  basis.   It  is,  at  the 
same  time,  increased  by  a  column  vector  corresponding  to  a  new  artifi- 
cial variable. 

The  new  (m+1)  x  (m+1)  basis  is: 

B    0 

By  partitioning  of  the  B  matrix,  the  new  inverse  is  readily 
obtai  ned: 


B 


1 


B 


1 


-Kb 


1 


1 


Upon  adding  a  new  constraint  to  a  problem,  one  of  three  cases  will 
result: 

(1)  the  new  constraint  is  satisfied; 

(2)  the  new  constraint  is  not  satisfied  and  the  value  of  the 
additional  artificial  variable  is  positive  in  the  basic  solution; 

(3)  the  new  constraint  is  not  satisfied  and  the  value  of  the 
additional  artificial  variable  is  negative  in  the  basic  solution. 

In  the  first  case,  the  new  constraint  has  no  effect  and  the 
optimal  solution  to  the  entire  original  problem  has  been  obtained  if 
no  further  constraints  are  to  be  added.  The  optimal  value  of  z  is  not 
affected  by  the  new  constraint. 

In  the  second  case,  the  two-phase  revised  simplex  procedure  is  used 
to  first  drive  out  the  artificial  variable  and  then  to  maximize  z. 
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In  the  third  case,  if  the  new  artificial  variable  is  assigned  a  zero 
cost  in  the  z  equation,  the  dual  simplex  procedure  can  then  be  applied  to 
the  infeasible  primal  to  obtain  an  optimal  solution  to  the  dual,  and 
hence  an  optimal  solution  to  the  primal.  The  artificial  variable  is  con- 
sidered to  be  a  legitimate  variable  of  the  original  problem  in  the  dual 
simplex  approach.  As  a  consequence,  it  never  appears  in  the  w  equation. 
For  convenience,  we  will  carry  along  the  "w"  row  of  the  tableau,  for  use 
when  later  constraints  are  added. 

When  a  solution  is  obtained  which  satisfies  the  added  constraint, 
the  optimal  solution  to  the  original  problem  has  been  found  if  there  are 
no  new  constraints  to  be  added. 

The  addition  of  a  new  constraint  to  a  linear  programming  will  have 
the  effect  of  either  decreasing  the  previously  obtained  maximum  solution 
or  leaving  it  unchanged.  That  is,  letting  the  subscript  on  z  denote  the 
number  of  constraints, 

max  zm+]  <  max  zm  . 

Solution  of  a  linear  programming  problem  by  step-by-step  addition 
of  constraints  may  result  in  one  or  more  unbounded  solutions  to  the  sub- 
problems  if  the  initial  constraints  have  fewer  variables  than  are   in- 
cluded in  the  objective  function.  This  causes  no  difficulty,  however, 
as  the  addition  of  one  or  more  new  constraints  will  serve  to  place 
bounds  on  the  problem,  unless  the  original  problem  is  unbounded. 

The  method  of  step-by-step  addition  of  constraints  has  the  advan- 
tage that  an  infeasible  solution  at  an  early  stage  will  determine  that 
the  original  problem  has  no  feasible  solution,  and  no  new  constraints 
need  be  added.  The  solution  procedure  is  then  terminated. 
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k.      Sample  Problem 

Consider,  as  an  example,  the  problem: 
maximize 

z  =  2x.  +  kx 
subject  to: 

3x1  +  kx~   +  6x3  <  2k 

ifx1  +  3x2  +  12x,  <  2k 

x,  +  x2  +  4x^  =  8 


and, 


x.  >  0  for  j  =  1,  2,  3. 


Adding  slack  variables  and  rewriting  the  problem,  we  have 
maximize 

subject  to: 

z  -  2x        -  kx  =  0 

3x1  +  4x2  +  6x0  +  x.      -2k 

4x]  +  3X2   +  12Xo        +  Xr    =    2k 

Xj  +  x2  +  kx-,  =   8 

and, 

x •  >  0   for  j  =  1 ,  ..  ,  5 

For  our  solution  by  the  SSAC  procedure,  the  modified  w  equation 
for  the  first  sub-problem  will  be: 

w  -  3xj  -  A-x2  -  6x0  -  x^  -  -2k 
because  the  original  w  equation  is 

w  +  Xs    -    0, 
where  X£  is  the  artificial  variable  introduced  into  "Rl". 
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The  initial  tableau  will  bet 

2   w   x6  x1   x2  x3  xk     x5  ^ 


lz.l 

1 

1 

1 

-2 

-3 

0  -k        0 
-i+  -6  -1 

0 
0 

0 
-2k 
2k 

'RT' 

3 

4  ©  1 

0 

W 

/K 

-e-  =  ik/G  =  4 


In  Phase  I  we  will  maximize  w.  We  choose  the  most  negative  value  in 
the  modified  w  equation  and  pivot  on  the  element  which  meets  the  mini- 
mum -©•  criterion  according  to  the  usual  (primal)  simplex  procedure. 
Note  that  in  the  first  tableau  the  first  pivot  will  always  be  such  as 
to  drive  out  the  first  artificial  variable. 

We  pivot  using  the  product  form  of  the  inverse.  1 5 1  The  ft     vector 
corresponding  to  the  Xo  column  will  be 


1 


k/6 

2/3 

6/6 

= 

1 

j/e 

1/6 

and, 


1 

0  2/3 

0 

1    1 

0 

0   1/6 

E  = 


_=1     -1 
The  new  inverse  is  determined  from  B   =  EB   as  follows 


_-l 

B 


1 

0  2/3 

0 

1      1 

0 

0  1/6 

1 

0 

0 

0 

1 

0 

= 

0 

0 

0 

1   0  2/3 
0   0    0 
0   0   1/6 
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The  new  tableau    is   now 

Z    W    X£   X,   X2   X~   X^   Xp   b 


1 

0   2/3 

0  8/3 

0  2/3 

0 

16 

0 
0 

1    1 
0  1/6 

— i 

0  0 

0  0 

0 

0 

k 

1 

— , 0 

where  x_  is  now  a  basic  variable;  all  values  of  z^-c.-  >  0,  and  we  have 
achieved  the  "first"  optimal  solution. 

To  add  a  new  constraint,  the  "original"  tableau  is  augmented  by 
one  row  and  the  column  vector  for  x-,,  where  Xj   is  the  artificial  variable 
associated  with  "R2".   Ignoring  the  "w"  row  for  the  time  being,  the 
augmented  "original"  tableau  will  be 

Z     W     X£    Xy   X^    X£    Xo    X^    Xr    b 


1 

-2 

0 

-k 

0 

0 

0 

1 

3 

k 

6 

1 

0 

2k 

1 

k 

3 

12 

0 

1 

2k 

We  compute  the  X   vector  by  determining  from  the  previous  solutions 
which  vectors  are  in  the  basis.-  Observe  that  the  coefficient  vectors 
associated  with  z  and  w  will  always  be  in  the  basis  since  we  are  maxi- 
mizing z  and  w.  The  first  two  elements  of  the  Q   vector  will  therefore 
always  be  zeros.  The  third  element  in  this  case  will  be  the  coefficient 
of  x   in  the  new  constraint. 


)T  =  [°   °    12] 
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The  augmented  inverse  of  the  basis  is  determined  by  partitioning: 


_-1 
B 


-tfV1      1 


The  product  -  Q   B~      ,  in  this  case,  is 

1    0  2/3 

To    0   -  1 2  1 


0    1    1 
0   0   1/6 


[o   o  -2]   , 


and  the  new  inverse  becomes 


1 

0 

2/3 

0 

_-1 

B 

0 
0 

1 

0 

1 
1/6 

0 
0 

0 

0 

-2 

1 

Our  next  step  is  to  determine  if  the  value  of  the  new  artificial  variable 
is  positive  or  negative  in  the  basic  solution.   We  find 

x7  =  2k   -  12x3  =  -2k. 

Because  it  is  negative,  our  procedure  tells  us  to  assign  Ci  =  0  to  the 
new  artificial  variable,  Xy,  and  use  the  dual  simplex  algorithm  to  drive 
out  the  artificial  variable. 

z   w   xg  x-j     x.   x»   x,  x.   Xj.   b 


1 

0 

0  2/3   0 

1  1    0 

0  8/3 

0  2/3  0 

0 

0 

0   1/6   0 

1/2 

0 

0   -2    1 

-* #-1 

©-5 

0  -2    1 

—4 , 1 

-2k 
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We  pivot  on  the  row  having  the  most  negative  b-  element  in  the  dual 
simplex  method.  As  a  consequence,  the  new  artificial  variable  will  be 
dropped  as  a  basic  variable. 

Our  pivot  element  is  determined  for  a^ •  <  0. 


■e-  = 


=  Min    ["-a^/a^j]   =   Min    lo,   8/15,    2/6 1  =   0, 


The   new   inverse    is   determined   by   the   product    form  of   the    inverse. 


1 


0 
0 

»A 

-1/2 


E    = 


10        0        0 
0         10        0 
0        0        11/4 
0        0        0-1/2 


The  current  tableau  becomes: 


z 

w        x6     x7 

x]      x2     x3     x4     xc 

b 

1 

0 

0  2/3      0 

1  1         0 

0     8/3      0     2/3      0 

0 

0 

0   -1/3    \/k 

0   -7/12    1(^/3)   ]/k 

-2 

0 

0        1-1/2 

0 ■ ■'■■■' 

1                ^       1 
U $ 

12 

-9-=  Min  [32/7,  2]  =  2 
Since  the  up-dated  requirements  vector,  b,  still  has  a  negative  com- 
ponent, the  dual  simplex  algorithm  must  again  be  applied.   Pivoting  on 
the  x^  column  gives; 

x6  x7   xl   x2   x"3  xk     x5   k 


w 


1  0   0   1/2 

0  110 

0  0   1  -3A 

0  0   0   1/4 

^ % 


0  3/2  2   0   1/2 


12 

0 
6 
6 
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Since  all  elements  of  b  are  now  non-negative,  our  solution  is  optimal 
Notice  that  the  "w!l  row  is  retained  but  not  operated  on  in  the  use  of 
the  dual  simplex  approach. 

Adding  the  third  constraint,  the  augmented  "original"  tableau 
becomes 

Z     W     X£   Xy    Xg   X^    X2    Xo    X^    Xr     b 


1 

-2 

0 

_z+ 

0 

0 

0 

1 

3 

k 

6 

1 

0 

2k 

1 

k 

3 

12 

0 

1 

2k     • 

1 

1 

1 

k 

0 

0 

8 

Again,  the  "w"  row  is  ignored  until  after  we  determine  if  the  new 
artificial  variable,  Xn,  is  positive  or  negative  in  the  basic  solution. 
After  bringing  in  the  new  constraint, 

#  =  [o   0   0   1  ] 
and, 

-^B"1  =  [o      0   0  -i/y    . 

The  value  of  the  new  artificial  variable,  xg,  for  constraint  "R3",  is 

x8  =  8  -  x,  =  2. 
Since  the  value  of  xo  is  positive  in  the  basic  solution,  we  can  use 

Phase  I  of  the  revised  simplex  procedure  to  move  to  the  optimal  solution 

for  the  complete  "original"  problem.   The  augmented  "w"  equation  in  the 

"original"  tableau  will  be 

z   w   x6  x-j     xg  x1   x2  x^  x^  x^   b 
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and  the  current  tableau  becomes 

Z     W     X/-   X-,   Xn    X,    X-    X,    X^   X ,- 


1 

0 

0      1/2      0 

2 

0 
0 

1 

0 

1         0        0 

i  -3A    o 

-1      -1       -if        0        0 

6 

-3 

0 

0 

0      \/k      0 

3 

6 

0 

-*- 

0 
— * — 

0-1/4      1 

&- 

© 

— # 1 

2 

«-=  Min  T6/3,  2/1 1  =  2 
Pivoting  on  the  appropriate  element  of  the  column  having  the  most  nega- 
tive value  in  the  "w"  row,  we  bring  x,  into  the  basis.  The  resulting 
tableau  is  then 


w 


:6  x7  x8  x1  x2  x3  xk     x5 


1 

0 

0        1      -2 

0 

0 
0 

1 

0 

1      -1        k 
1    -3/2      3 

-1         0 

0 

0      -1 

0 

0 

0 

0        1      -3 

Q 

0 

0 

' — •- 

0 

— • — 

0  -\fk      1 

0 

#~ 

-#» 

Since  all  of  the  coefficients  of  the  "w"  row  are  not  non-negative, 
we  again  pivot  on  the  column  having  the  most  negative  element,  or  in 
the  case  of  ties,  the  left-handed  one  of  the  tied  columns.   Pivoting 
on  the  Xi  column,  we  bring  Xi  into  the  basis. 

After  this  iteration,  we  find  that  all  elements  of  the  "w"  row  are 
zero  and  at  the  same  time,  all  z-  -  c-  >  0  so  we  have  completed  both 
Phase  I  and  Phase  II.  The  final  tableau  is: 
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z 

w 

x6     x7 

x8 

x1 

x2 

x3 

xk 

x5 

b 

1 

0 

0        1 

-2 

0 

1 

0 

0 

1 

8 

0 
0 

1 
0 

1         0 
1   -3/2 

1 

3 

0 

0 

0 

0 

0 

0 
12 

1 

0 

0 

0        1 

-3 

1 

0 

0 

0 

o  -iA 

1 

1 

2 

We  have  now  arrived  at  the  optimal  solution  of  the  "original"  problem: 

z  -  8 
x,  =  0 
Xo  =  2 
xk   =    12 


and  all  other  x-  =  0. 


2k 


5.   Programming  Technique 

As  a  means  of  testing  the  feasibility  of  solving  general  linear 
programming  problems  by  the  SSAC  method,  the  solution  technique  was 
coded  in  FORTRAN  IVg  for  use  on  the  IBM  360/67  computer. 

One  subroutine  was  designed  which  would  carry  out  the  solution 
procedure  by  either  the  two-phase  revised  simplex  method  or  the  dual 
simplex  method  as  appropriate.  This  same  subroutine  was  used  to  solve 
the  problems  by  both  the  standard  revised  simplex  procedure  and  by  the 
SSAC  procedure.   By  using  the  same  subroutine  for  both  methods,  it  was 
hoped  that  any  bias  which  might  result  from  programming  technique 
could  be  avoided.  A  driving  routine  was  designed  which  would  first 
solve  a  problem  by  the  revised  simplex  method  and  then  re-solve  the 
same  problem  using  the  SSAC  method.  The  two  sections  of  the  program 
were  then  timed.   Sections  of  the  program  not  germane  to  the  method 
being  investigated,  such  as  the  reading  in  and  printing  out  of  data, 
were  not  included  in  the  timing.  The  number  of  iterations  required 
for  solution  by  each  method  was  tabulated. 


'The  timing  routine  was  developed  by  Lt  E.  A.  Singer,  a  student 
at  the  Naval  Postgraduate  School. 
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6.   Efficiency  of  the  Algorithm 

Several  small  problems,  for  which  hand  solutions  could  easily  be 
obtained,  were  used  for  preliminary  testing  and  debugging  of  the  program, 
A  number  of  larger  problems  were  then  solved  to  obtain  a  limited  ex- 
perimental verification  of  the  new  procedure.  The  problems  were  chosen 
from  three  categories  of  problem  types;  mixing  problems,  transportation 
problems,  (including  a  network  problem  and  a  transshipment  problem),  and 
caterer  problems.   Results  based  on  this  preliminary  comparison,  as 
shown  in  Table  I,  were  inconclusive. 

In  order  to  obtain  the  solutions  to  a  large  number  of  problems,  and 
as  a  means  of  avoiding  the  considerable  time  and  effort  required  to 
input  data  by  hand,  a  routine  was  designed  which  would  generate  random 
problems.  This  routine  utilized  a  random-number  generator  to  generate 
elements  for  the  A  and  b  matrices.  To  insure  the  existence  of  a 
bounded  optimal  feasible  solution,  the  problems  were  formulated  as; 
Maximize 


subject  to: 


z  =  ex 


Ax  >  b 


and, 

Cj  <0 

x .  >  0      for  j  =  1 ,  . . .,  n 

ai  >° 

b.  >  0      for  i  =  1 ,  . . .,  m 

The  problems  were  generated  to  have  70  variables,  including  slack 

variables,  and  20  constraints.  The  distribution  of  the  coefficients 

was  uniform  over  the  following  intervals; 
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a . .  ,   uniform  (0.1); 

b.   ,   uniform  (0,5); 

c-   ,   uniform  (-1,0). 
A  total  of  forty-six  problems  were  solved  using  this  method  of  problem 
generation.   A  tabulation  of  solution  times  for  these  random  problems 
is  given  in  Appendix  I. 

A  comparison  of  the  solution  times  of  these  forty-six  problems 
shows  that  the  method  of  step-by-step  addition  of  constraints  was  faster 
in  thirty-four  cases.   The  mean  solution  time  for  the  SSAC  procedure 
was  3.9^  seconds  faster  than  the  mean  solution  time  by  the  revised 
simplex  procedure.   By  applying  an  appropriate  statistical  test  to  the 
solution  results,  it  was  determined  that,  with  95%  confidence,  the  mean 

difference  in  solution  times  for  the  two  methods  is  not  less  than  2.27 

2 
seconds.   Therefore,  it  can  be  concluded  that  the  method  of  SSAC  is 

significantly  faster  than  the  revised  simplex  method. 


Detailed  computations  for  the  t-test  and  computation  of  lower 
confidence  limits  on  the  mean  solution  time  difference  are   given  in 
Appendix  I. 
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TABLE  I 


SOLUTION  TIME  RESULTS  OF  PROBLEMS  USED  IN  PRELIMINARY  INVESTIGATION 


Problem 

Revised  Simplex 
(RS) 

Add.  of  Constr. 
(SSAC) 

T  i  me  Ratio 

Iteratio 

ns  Time 
(sec) 

i 
Iterati o 

ns  Time 
(sec. ) 

SSAC/  R.S. 

fixing  Problems 

9 

.106 

15 

.188 

1.17 

laugh's  Diet 

Gasol  i  ne  BlencT 

10 

.766 

29 

l.i+46 

1.88 

Transportation  P 

r obi  ems 

118.829 

89 

51.352 

.43 

3  by  33  tableau 

85 

7  by  7  tableau 

41 

6.398 

38 

4.071 

.64 

3  by  5  tableau 

18 

.852 

22 

.700 

.82 

3  by  4  tableau 

13 

.446 

21 

.479 

1.07 

Transshi  pmentr 

18 

4.223 

20 

3.117 

.74 

Network  F low 

17 

3.318 

24 

2.507 

.75 

Caterer  Problems 

[1] 

Wardroom  Napkin 

35 

10.768 

39 

6.177 

.57 

[5] 

Hadley  Napkin 

10 

.745 

25 

.918 

1.23 
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7.   Concluding  Remarks 

It  is  important  to  emphasize  that  the  solution  time  is  a  function 
of  the  way  the  computer  program  is  written.   In  the  code  employed,  a 
full  set  of  artificial  variables  was  generated  for  each  problem.  An 
algorithm  which  takes  advantage  of  existing  slack  variables  for  an 
initial  feasible  solution  might  well  prove  to  be  faster  than  the 
present  program. 

Round-off  error  and  exponent  underflow  can  greatly  affect  the 
solution  technique.  The  use  of  the  product  form  of  the  inverse  in  the 
pivoting  operation  relieves  this  situation  somewhat.   In  the  experi- 
mental algorithm  a  routine  was  employed  which  set  any  element  having 
an  absolute  value  less  than  .0001  equal  to  zero.  A  better  method  might 
be  to  use  double  precision  mode  for  computations  and  then  allow  values 
smaller  than  .0001  to  be  carried  along  in  the  solution. 

By  using  a  single  subroutine  for  the  simplex  iteration  procedure 
in  both  methods  of  problem  solution,  it  was  hoped  that  any  inconsistency 
due  to  programming  technique  could  be  kept  to  a  minimum.  That  is, 
necessary  computations  for  the  iteration  procedure  were  carried  out  in 
the  same  sequence  for  both  methods  of  problem  solution. 

It  is  recognized  that  the  problems  selected  are  not  necessarily 
a  representative  sampling  of  linear  programming  problems.  The  manual 
input  problems  were  selected  primarily  because  they  were  large  enough 
to  allow  the  step-by-step  addition  of  constraints  to  be  demonstrated, 
yet  small  enough  to  handle  conveniently  as  data  inputs.   The  randomly 
generated  problems  were  considered  to  be  of  a  size  which  was  large 
enough  to  effectively  test  the  procedure,  subject  to  available  com- 
puter time. 
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The  results  of  this  limited  number  of  tests  indicate  that  the  SSAC 
method  is  worthy  of  further  investigation.  This  further  effort  could  be 
directed  in  one  of  several  areas.   The  present  method  should  be  applied 
to  a  large  number  of  more  diverse  problems  in  order  to  obtain  a  better 
data  base  for  verification  of  the  results  already  obtained,  and  to  deter- 
mine more  accurately  the  advantage  of  this  method  over  the  revised 
simplex  procedure.   At  the  same  time,  it  would  be  possible  to  determine 
some  bounds  of  effectiveness  of  this  procedure  as  to  the  size  and  struc- 
ture of  problems. 

Modification  of  the  step-by-step  procedure  might  be  attempted  to 
take  advantage  of  an  existing  basis  in  the  original  problem  so  as  to 
require  generation  of  fewer  artificial  variables.  An  attempt  might 
also  be  made  to  modify  the  step-by-step  addition  of  constraints  pro- 
cedure for  application  to  the  primal-dual  algorithm. 
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APPENDIX  I. 

Statistical  Testing  of  Random  Problem  Solutions 

If  it  can  be  assumed  that  the  solution  times  obtained  by  each 
method  form  two  normal  distributions,  we  can  test  to  determine  if  the 
mean  solution  time  by  the  addition  of  constraints  method,  /kSAC  1#  s 
less  than  the  mean  solution  time  by  the  revised  simplex  method,   /-"rc 

Letting  X.  be  the  solution  time  obtained  by  the  revised  simplex 
procedure,  and  Y-  the  solution  time  obtained  by  the  addition  of  con- 
straints procedure  for  problem  i,  a  pairwise  comparison  of  the  results 
of  the  two  methods  can  be  made. 

Consider  the  hypothesis  that  the  mean  solution  time  difference, 

/"*  ~  /''RS  "  /"SSAC*  1S  non-positiv6*  That  is, 

H0:   yu  <  0, 

with  the  alternative  hypothesis, 

Hj:   yu  >  0  . 
To  show  that  the  method  of  addition  of  constraints  is  faster,  we  must 
be  able  to  reject  Hq. 

If  we  form  a  t  statistic  for  n  =  k6   samples, 


t  =  D  y   n 
S 


where, 


D  =_L  X  (X.  -  Y.) 

n   i=1 
and, 

S2  =   1     Z  (D.  -  D)2   , 

n-1    i=1    ' 

and  then  compare  it  against  tabulated  values  of  the  cumulative  t  dis- 
tribution, we  would  reject  the  hypothesis  if 


4  £  V"-i  • 
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In   this    case  we    choose   a    100(1    -0( )    =   95%   confidence    level.      The   com- 
putation   of   the   t    statistic    for    n    =  46,    S    =   6.74   gives, 


t   -     3.94     7  46  =.  3.97      . 
6775 

The  tabulated  value  for  %,„    __,  ,._  s  1.68.   Since  t  V  t  ,    .we 

(0.05), 45  '        (^,n-l 

reject  the  hypothesis  that  the  revised  simplex  method  is  faster. 

A  lower  confidence  limit  on  the  mean  difference  in  solution  times 
can  be  determined  from  the  expression 
Pre 


ob    |   (D    -  /jl)  Yn~  ^    t^l    =    1    -  0<  . 

confidence 

(d  -  >iL)yr  <  tfl( 


Since  with    95%  confidence 


we   can   solve   for   the    lower    limit    on     jj, 

For    the    given  data, 

3.94   -    (1.68)(6.74)     </  yiL 

(6.785) 

or, 

2.27  <yjL  . 

So  we  can  say  with  95%  confidence  that  the  mean  difference  in  solution 
times  is  no  less  than  2.27  seconds.  That  is,  the  mean  solution  time 
obtained  by  the  step-by-step  addition  of  constraints  is  at  least  2.27 
seconds  less  than  the  solution  time  obtained  by  the  revised  simplex 
procedure. 
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TABLE 

II 

SOLUTION 

TIME 

RESULTS    OF    RANDOMLY 

GENERATED    PROBLEMS 

Prob. 

RS 

SSAC 

No. 

Time    (X^ 
(sec. ) 

Iter. 

Time   (Y-) 
(sec. ) 

Iter . 

X-    -   Y. 

(Xi    -  Yi)2 

1 

23.90 

68 

5.11 

25 

18.79 

220.5 

2 

18.19 

51 

7.08 

56 

11.11 

51.4 

3 

14.97 

42 

IO.39 

71 

4.58 

0.4 

4 

28.25 

82 

18.34 

101 

9.91 

35.6 

5 

20.21 

58 

16.66 

92 

3.55 

0.2 

6 

23.31 

64 

9.88 

64 

13.43 

90.1 

7 

27.24 

79 

25.23 

179 

2.01 

3.7 

8 

12.10 

33 

5.14 

24 

6.96 

9.1 

9 

13.47 

37 

17.63 

93 

-4.16 

65.6 

10 

13.47 

37 

5.56 

38 

7.91 

15.8 

11 

13.13 

36 

18.50 

113 

-5.02 

80.3 

12 

14.50 

40 

5.69 

38 

8.81 

23.7 

13 

9.72 

26 

14.68 

80 

-4.96 

79.2 

14 

13.48 

37 

7.68 

59 

5,80 

3.5 

15 

13.14 

36 

10.10 

91 

3.04 

0.8 

16 

10.76 

29 

18.43 

93 

-7.67 

134.8 

17 

15.53 

43 

15.41 

90 

0.12 

14.6 

18 

19.62 

56 

7.29 

48 

12.33 

70.4 

19 

24.38 

69 

12.05 

65 

12.33 

70.4 

20 

19.16 

54 

10.30 

51 

8.86 

24.2 

21 

22.23 

63 

8.68 

56 

13.55 

92.4 

22 

21.93 

62 

10.83 

60 

11.10 

51.3 

23 

26.35 

75 

19.69 

87 

6.66 

7.4 

24 

23.65 

67 

7.69 

45 

15.96 

144.5 

1  25 

14.40 

43 

5.22 

38 

9.18 

27.5 

26 

9.72 

28 

12.20 

90 

-2.48 

41.2 

27 

14.73 

44 

7.02 

52 

7.71 

14.2 

28 

10.34 

30 

10.93 

70 

-0.59 

20.5 

29 

11.58 

34 

14.32 

97 

-2.74 

44.6 

30 

11.62 

34 

9.06 

58 

2.56 

1.9 

31 

12.51 

35 

11.81 

81 

0.70 

10.5 

32 

12.17 

34 

5.43 

32 

6.74 

7.8 

33 

13.16 

37 

5.44 

38 

7.72 

14.3 

34 

11.83 

33 

5.55 

40 

6.28 

5.5 

35 

12.18 

34 

17.25 

93 

-5.07 

81.2 

36 

12.51 

35 

16.03 

73 

-3.52 

55.7 

37 

11.18 

31 

14.30 

101 

-3.12 

49.8 

I  38 

15.18 

43 

15.07 

90 

0.11 

14.7 

139 

9.24 

27 

7.90 

54 

1.34 

6.8 

Uo 

11.98 

36 

12.72 

81 

-0.74 

21.9 

41 

11.06 

33 

9.99 

74 

1.07 

8.2 

42 

11.41 

34 

7.86 

62 

3.55 

0.2 

k3 

13.52 

41 

6.88 

57 

6.64 

7.3 

kk 

11.36 

34 

8.63 

51 

2.73 

1.5 

45 

11.68 

35 

24.94 

142 

-13.26 

295.8 

k6 

9.53 

28 

8.16 

58 

1.37 

6.6 
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APPENDIX  II. 


Flow  Diagrams  of  the  Computer  Progra 


m 


No 


C  MAIN     J 


X 


.Read  in 
Data 


Store  Input  in 
Temporary  Arrays 
for  Computations 


Compute  "w"  Row 


Call  SIMPLX 


PRINT  Results  from 
Revised  Simplex  Method 


-t 


Reset  Problem  for 

Solution  by 

Addition  of  Constraints 


I 


Add  New  Constraint 


Compute  0 


Compute  New  Inverse 


■& 


Call  SIMPLX 


-t 


Have  All  Constraints 
Been  Added? 


Yes 


PRINT  Results  from 
Addition  of  Constraints 
Method 

(     END     ) 
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2.   Simplex  and  Dual  Simplex  Iteration  Subprograr 


C  Sub  routine  "N 
SIMPLX    J 

\  _ 


Compute  Elements 
of  "w"  row 


Is  New  Art  if icia 1 
Negative  in  Basis? 


No 


Is  Phase  I  Complete? 


No 


Find  Min  Negative 
Element  of  "w"  Row 


Find  Minimum  -$ 


I 


Compute  h  Vector 


Compute  New  Inverse 


Is  Phase  I  Complete? 


Yes 


X 


Is  Phase  II  Complete? 


Yes 


Yes 


-&-I  1000 


Yes 


c/"  RETURN  J 


No 


Find  Min  Negative 
Element  of  "z"  Row 
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(Dual  Simplex  Section) 


Find  Most  Negative 
Element  of  b  Vector 


m 


nd  M  i  n  i  mum  •$■ 


Compute  ^1  Vector 

Z3 


Compute  New  Inverse 


£ 


Is  Phase  II  Complete? 
ii  " 


c 


Yes 
! 


RETURN 


D 


No 
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